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ABSTRACT 

The role of MHD turbulence in the cosmic ray acceleration process in a volume with 
a reconnecting magnetic field is studied by means of Monte Carlo simulations. We 
performed modelling of proton acceleration with the 3D analytic model of station- 
ary reconnection of Craig et al. (1995) providing the unperturbed background condi- 
tions. Perturbations of particle trajectories due to a turbulent magnetic field compo- 
nent were simulated using small-amplitude pitch-angle momentum scattering, enabling 
modelling of both small and large amplitude turbulence in a wide wave vector range. 
Within the approach, no second-order Fermi acceleration process is allowed. Compar- 
ison of the acceleration process in models involving particle trajectory perturbations 
to the unperturbed one reveals that the turbulence can substantially increase the 
acceleration efficiency, enabling much higher final particle energies and flat particle 
spectra. 

Key words: magnetic field: reconnection - cosmic rays - acceleration of particles - 
MHD Sun: flares - turbulence 
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1 INTRODUCTION 



X 



Particle acceleration processes in regions of magnetic field 
reconnection can play an important role in various astro- 
physical sites by providing high energy particles, heating 
thermal plasma and by changing the magnetic field con- 
figuration (e.g. Somov 1994). The first of these processes, 
accelerating particles to cosmic ray energies, acts due to 
electric fields occurring in the reconnection region. Two ap- 
proaches were applied to describe the acceleration process. 
In the more basic one, analytic or semi-analytic models for 
the reconnection region provided background for the deriva- 
tion of particle trajectories and discussion of the energy dis- 
tribution of particles escaping from the reconnection vol- 
ume. A simple discussion of particle trajectories presented 
by Speiser (1965; see also Sonnerup 1971) shows that even 
a tiny magnetic field component near the neutral current 
sheet efficiently removes ergetic particles and decreases its 
final energy. Further studies by Stern (1979) and Wagner et 
al. (1981) revealed the importance of the O-type neutral line 
regions for the acceleration, due to partial particle trapping 
in the accelerating volume. Deeg et al. (1991) and Somov & 
Kosugi (1997) discussed several features of the acceleration 
process in the X-type stationary reconnection. The time de- 
pendent effects of the acceleration process were introduced 
by considering the background magnetic and electric fields 
as derived from simple MHD simulations (cf. Sato et al. 



1982, Scholer & Jamitzky 1987, Atkinson et al. 1989, Ze- 
lenyi et al. 1990). 

Substantial progress in considering realistic reconnec- 
tion processes came from the two dimensional (= 2D; 
Matthaeus et al. 1984, Ambrosiano et al. 1988, Scholer & 
Jamitzky 1989, Veltri at al. 1998) and 3D (Birn & Hesse 
1994, Schopper et al. 1999, Kliem et al. 1998) MHD mod- 
elling involving a perturbed magnetic field as the initial 
condition. This approach yields a complex structure with 
X-type and O-type null points and magnetic field pertur- 
bations of scales comparable to the macroscopic structure 
dimensions. Ambrosiano et al. (1988) show that such a tur- 
bulent neutral point mechanism influences the acceleration 
process in two ways. It enhances the reconnection magnetic 
field while producing a stochastic electric field that gives 
rise to momentum diffusion, and it also produces magnetic 
bubbles and other irregularities that can temporarily trap 
test particles in the strong reconnection field for times com- 
parable to the magnetofluid characteristic time. As a result, 
the very flat particle spectra formed can extend to higher 
energies in comparison to the unperturbed conditions. 

One should also note that a realistic description of the 
acceleration process requires considering the full 3D config- 
uration of the reconnection region (e.g. Birn & Hesse 1994, 
Schopper et al. 1999, Veltri et al. 1998). Besides the role 
which the mean field structure plays, particle transport de- 
pends in a qualitative way on the turbulence dimensionality 
(cf. Giacalone & Jokipii 1994, Michalek & Ostrowski 1996) 
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and the involved wave modes (e.g. Schlickeiser & Miller 1997, 
Michalek & Ostrowski 1998). 

Besides providing a substantial improvement in under- 
standing of the process, a general deficiency of the above 
mentioned papers discussing particle acceleration in the tur- 
bulent reconnection regions are the complicated interrela- 
tions between such factors as the assumed or derived mean 
magnetic field structure, its dimensionality, the assumed or 
derived form of turbulence, and the particle acceleration pro- 
cess. For example in the paper of Ambrosiano et al. (1988) it 
was difficult to separate the role of diffusive particle motions 
from trapping by the structures formed in/near the recon- 
nection layer. In order to clarify the situation, with only one 
of these factors modifying the acceleration process in turbu- 
lent reconnection, we decided to use a simple model which 
allows us to consider the role of random particle trajectory 
perturbations separately. The process acting in the vicin- 
ity of the X-type reconnection does not provide a means for 
particle trapping in the mean (unperturbed) structure of the 
reconnecting field. As described in the next section, we use 
the 3D analytic model of stationary reconnection of Craig 
et al. (1995) to define the background model, being a refer- 
ence for the model involving particle pitch-angle momentum 
scattering due to MHD turbulence. Perturbations of parti- 
cle trajectories due to a turbulent magnetic field compo- 
nent were simulated using small-amplitude pitch-angle mo- 
mentum scattering, enabling modelling of both small and 
large amplitude turbulence in a wide wave vector range. No 
second-order Fermi acceleration processes are allowed within 
this approach. In section 3 we present results of particle 
energy spectrum modelling for different scattering ampli- 
tudes. Comparison of the acceleration process to the unper- 
turbed one confirms that the turbulence can substantially in- 
crease the acceleration efficiency, enabling particles to form 
flat high-energy spectra with much increased final energies. 
Then, in section 4, we briefly summarize these results. 

In the present simulations we scale the parameters of 
the reconnection region to those characteristic for the solar 
flare (e.g. Miller et al. 1997). However, we do not aspire to 
present a flare acceleration model. We provide all numerical 
values in SI units. 



2 SIMULATIONS OF THE ACCELERATION 
PROCESS 

2.1 The 3D reconnecting magnetic field structure 

The considered reconnection region (Fig. 1) involves in- 
flowing plasma from the 'top' and the 'bottom' toward the 
x — plane, and outflowing to the sides. All physical quan- 
tities, the flow velocity, the magnetic field and the electric 
field are symmetric with respect to the origin of the refer- 
ence frame, with the reconnection "plane" (magnetic field is 
strictly only in a single point in the center) being slightly 
inclined with respect to the x = plane. To describe this 
structure analytically we use a 3D model of the stationary 
reconnection by Craig et al. (1995), which provides the an- 
alytic structure for the magnetic field B, the electric field 
E and the plasma velocity V . The following parameters are 
introduced in this model: a, A, r\ and k define topologic 
structure of the magnetic field; p define the reconnection 
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Figure 1. Upper panel: Magnetic field vectors in three layers 
y = — 3e7, y = 0, y = 3e7 of the considered recconection volume. 
The y component of magnetic field vector has been increased 6 
times to demonstrate of the stricture. Middle and bottom panels: 
The magnetic field structure (full lines) and the velocity field 
(dashed lines) in the y = plane of the considered reconnection 
region. In both these panels the border of the figure also represents 
the particle escape boundary. One should note different scales 
along the axes in these panels. The small black rectangle in the 
middle shows the size of the region where particle scattering due 
to turbulence is imposed. In the lower panel, the structure of the 
unperturbed magnetic field and the velocity field in this small 
turbulent volume is presented. 
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area thickness; p is the plasma density; £ denotes the exter- 
nal electric field component; po is the magnetic permeability 
of the vacuum. With the use of these parameters and the no- 
tation involving a radius vector r — {x , y , z} and the unit 
vectors x, y and z along the respective axes, 



B = \P(f)+Q(x) 



(2.1) 



E = 



where 



V = \P(r) + \Q(x)] -L= , (2.2) 
y + £ [1 — 2xfi da,w(px)] z , 
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and 




daw(p:r)y + 



2/' 



Z'(0)erf(/xa;) + Z(0) 



(2.3) 



(2.4) 



(2.5) 



daw(x) = / exp(i 2 — x 2 )dt , erf(is) = / exp(— t 2 )dt 



p 2 = — 1-A . 
M 2r, 

The function Z(a;) of Craig et al. (1995) satisfies the differ- 
ential equation 

a{l- \ 2 )[(1~ k)Z + xZ']+ V Z" =aj 2 \(2-k)x , (2.6) 

where one requires 72A: = (72 or fc must be zero). With 
these expressions one can construct the reconnection region 
in two ways. It is possible to choose the parameters £ , Z'(0), 
Z(0), a, P, rj, fi in such a way that either the magnetic field 
component Bg or By vanishes in some area: 

• i.) By vanishes where Xaky = —7^7 da.w(jj,x) 

• ii.) B z vanishes where /3(1 - k)z = &Z'(0)eri(fix). 

The maximum value of r allowing for one of the above condi- 
tions to be satisfied gives the half length of the reconnection 
region L (a scale for the layer with negligible magnetic field 
near the null point). For example, for the case (ii) one has 



7rZ'(0) 
4/3(1 - fc)/i 



(2.7) 



In the case (i), where the function daw determines the be- 
haviour of the magnetic field, the gradient of the magnetic 
field is greater than in the case (ii). As a result the simu- 
lation time for the required accuracy is longer. Moreover, 
within the turbulent volume, at a large distance from the 
centre the function daw grows slower resulting in less effi- 
cient acceleration. For these reasons in the present simula- 
tions we consider solutions determined with the condition 
(ii) satisfied. We take the following values of the model pa- 
rameters: a = 1.0 • 10~ 9 , A = -, where /3 = 1.0 ■ 10 -10 , 



f] = 1.0 • 10" 6 , k = 0.15 , p = M p • 10 16 (M p proton mass), 
£ = 0, Z(0) = 0, Z'(0) = 7- 10" 4 . Values of these parameters 
were derived from fitting the model magnetic, electric and 
velocity field to the mentioned solar flare conditions, with 
the plasma density equal to 10 16 m -3 . 

Let us note that the simple analytic form of the con- 
sidered solution limits the range of physically acceptable 
parameters to provide the possible boundary values of V, 
B, E . In the present simulations we choose such bound- 
ary values for the solution at the edge of our considered 
volume to provide physical parameters near the reconnec- 
tion layer that are close to the ones estimated for the solar 
flares. In the vicinity of the reconnection area, just outside 
of the reconnecting current layer, the magnetic field induc- 
tion is assumed to be B = 1.68 • 10 -3 T and the respective 
gyroradius (gyroperiod) for a 1 MeV proton is r g = 87 m 
(T g = 4 • 10" 5 s), and r g = 3.4 km (T s = 8.1 • 10" 5 s) 
for a 1 GeV proton. The assumed linear dimensions of the 
full considered reconnection region are Ax = 2 ■ 10 4 km, 
Ay = 2 ■ 10 5 km and Az = 1.5 • 10 5 km. The size of 
the volume of the perturbed magnetic field is in the cen- 
tre Ax = 2- 10 3 km, Ay = 6 ■ 10 3 km and Az = 2 ■ 10 3 km. 
One should note that the selected particular conditions are 
not essential for the presented considerations. 



2.2 Monte Carlo modelling procedure 

In the simulations we used a Monte Carlo approach includ- 
ing a trajectory splitting technique to improve statistics at 
larger energies (cf. Ostrowski 1991). To derive particle spec- 
tra we recorded 'weights' of particles escaping from the sim- 
ulation volume in a given energy range. For every escap- 
ing particle a randomly^] chosen particle still active in the 
simulations was split into two identical particles, each with 
half the weight of the original particle. Then the trajecto- 
ries of both particles were followed, but due to the applied 
random momentum scattering they evolve in different ways. 
The maximum simulation time t ma x (tmax is the upper limit 
in Fig. 4) was chosen to be sufficiently large to be unreach- 
able by high weight particles. Thus any further increasing 
of tmax does not influence the resulting spectra in a visible 
way. 

In the simulations, we inject test particles with initial 
energy Eo = 2 MeV in the vicinity of the reconnection null 
point and we follow their trajectories by numerical integra- 
tion of particle equations of motion. The integration is com- 
pleted when a particle either crosses the boundary of the 
considered reconnection region ('escape boundary'), or the 
time limit tmax = 100 s is reached. Particles were scattered 
in a small region around the reconnection null point only. 
This region is shown in the two upper panels of Fig. 1 as 
black rectangles, and in expanded form in the lower panel. 
The trajectory computation times are much longer for low 
energy particles and this is the main reason why we start 
with initial proton energies substantially larger than the 
thermal energy; for the discussion of the injection problem of 
energetic solar flare particles one can consult the discussion 
in Miller et al. (1997). As the long integrations performed 



* With probability of selecting a given particle being proportional 
to its weight. 
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here require high accuracy we use a variable step fourth- 
order Runge-Kutta integration with the parameters chosen 
in such a way that any further accuracy increase does not af- 
fect the simulated trajectories. In the simulation we use only 
100 particles because the required high accuracy of integra- 
tion leads to extensive integration times. Let us remember 
that with the applied trajectory splitting technique the num- 
ber of particles forming the spectrum at different energies is 
the same. 

We performed a number of numerical tests of the code 
applied in the simulations. We checked by hand the deriva- 
tion of the resulting particle phase space co-ordinates in a 
few randomly chosen individual integration steps of the code 
algorithm. Then we derived particle trajectories within a few 
simple uniform magnetic fields oriented randomly with re- 
spect to the chosen reference frame. The results coincided 
within the numerical accuracy with the derived analytic tra- 
jectories. Finally, for the actual considered magnetic field 
structure in the reconnection volume we positively checked 
conservation of particle energy with the plasma velocity set 
equal to zero. 

In the simulations particles gain energy mostly in the 
vicinity of the reconnection layer when drifting in the V AB 
electric fields. Away from this region, while moving across 
the magnetic field gradient, particles can gain and lose en- 
ergy, but the mean energy change is small. Somov & Kosugi 
(1997) estimated the accelerated particle energy as a prod- 
uct of the magnetic field, the plasma velocity and the recon- 
nection area length. The energy gains derived in our unper- 
turbed reconnection model are in agreement with this esti- 
mate if we take the distance traversed by a particle within 
the reconnection area as the required length. The effect of 
particle escape from the reconnection volume has been dis- 
cussed by Speiser (1965) for a highly simplified reconnection 
model. We confirm his results showing that a small vertical 
magnetic field component within the reconnection layer in- 
creases the particle escape substantially. 

In the simulations we consider the particle scattering 
process only in a small volume containing the reconnec- 
tion null point, where the magnetic energy is dissipated (see 
Fig. 1). The real reconnection regions are expected to show 
analogous structures, with the turbulence amplitude grow- 
ing toward the centre, including null points of the magnetic 
field. This choice also resulted from the fact that the Craig et 
al. (1995) three dimensional model of the stationary recon- 
nection is not too realistic at large distances from the null 
point, as, for example, it involves the inflow velocity grow- 
ing without a limit when increasing distance from the centre. 
For the perturbed trajectories model we use a simple pitch- 
angle scattering approach (e.g. Ostrowski 1991) intended to 
model the particle scattering at MHD waves. In this case 
integration of the particle equations of motion is performed 
in the electromagnetic field defined by the unperturbed ana- 
lytic model, but trajectory perturbations are introduced ev- 
ery constant time interval At, when the particle momentum 
vector is randomly scattered within a narrow cone along its 
original direction. The scattering is performed in the local 
plasma rest frame and it conserves particle energy in this 
frame. In the present simulations we consider the uniform 
momentum scattering within a cone with the half opening 
angle equal to 11°. The perturbation intensity is controlled 
by changing At and it is characterized with H = k±/k\\, the 



ratio of the cross-field diffusion coefficient to the diffusion 
coefficient along the magnetic field. The value of H was de- 
termined in auxiliary simulations performed in the uniform 
magnetic field with the value characteristic for a region close 
outside the reconnection current sheet (see section 2.3). One 
should also note that decrease of the magnetic field toward 
the reconnection site leads to increasing the effective turbu- 
lence amplitude (i.e. the value of H; in the limit of B = 
we have H = 1.0), but in the present simulations the parti- 
cle gyroradius is always larger than the reconnection layer 
thickness near the considered X-type null point. 

As the mean scattering time At is assumed to be con- 
stant within a given simulation run, particle trajectories are 
affected by perturbations with intensity depending on par- 
ticle energy and the background magnetic field. For non- 
relativistic particles with constant angular velocities of their 
gyration movements, assuming constant At is equivalent 
to introducing scattering acts at constant gyrophase steps. 
Thus the resulting value of H does not depend on energy, 
as expected for the flat wave power spectrum F(k) oc fc _1 . 
This slightly unrealistic wave spectrum allows, on the other 
hand, evaluation of the role of diffusive effects for the same 
scattering amplitude at all considered particle energies. A 
discussion of a more realistic Kolmogorov wave spectrum 
within the finite wave vector range will be presented in the 
next paper (in preparation) . However, as such a wave spec- 
trum carries more energy in long waves, the results are ex- 
pected to show a transition from our low H results to the 
large H ones. 

In attempting to compare our simplified scattering 
model with the real turbulence with an amplitude of waves 
resonant for a particle of a given energy, 8B r (SB 2 /8-k w 
F(k r ) ■ k r , where k r = 2ir/r g ), one can refer to a qual- 
itative discussion of energetic particle diffusion presented 
by Drury (1983). With his scaling Ky oc (5B r /B)~ 2 and 
ki oc (8B r /B) 2 , the amplitude for resonance waves can be 
evaluated as SB r /B « K 1 ^ 4 . 



2.3 Derivation of the N = k±/k« parameter 

The respective values of K = k±/k\\ were derived in auxil- 
iary simulations involving the spatially uniform background 
magnetic field with the induction of 1.5- 10~ 3 T and particles 
with energies equal to the initial energy Eo = 2 MeV. Tra- 
jectories of a large number of particles were followed with the 
imposed scattering process involving the momentum angular 
scattering, uniform within a cone of half opening angle equal 
to 11° and with the cone axis directed along the original mo- 
mentum vector. The only parameter varying between the 
simulations was the time interval between successive scat- 
tering events, At. The resulting diffusion coefficients were 
derived from growing particle dispersions along the back- 
ground field (ft||) and along two orthogonal axes perpendic- 
ular to the background field (= along the 1- or 2-axis). The 
results of such computations are presented in Fig. 2. Presen- 
tation of two derived values of ni/fty and K2/k\\ allows one 
to evaluate the accuracy of these computations. The values 
used in the paper are fits to the asymptotic (T — > oo) value 
of K = 0.5 (ki + K2)/k\\ ■ For a sequence of scattering times 
At = 10~ 6 , 10~ 5 , 10~ 4 and 10~ 3 we derived the respective 
values of N = 1.2- 10~ 2 , 6-10" 3 , 6-10" 5 , and6-10~ 7 . The re- 
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Figure 2. The diffusion parameter X = k,±/k,\ calculated during 
particle distribution evolution. T is the simulation time. For given 
initial parameters one has two resulting values of N for k± = K\ 
and k± = K2, presented as separate lines. 
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Figure 3. The resulting energy spectra for protons injected at 
E = 2 MeV: Model A: H = (the unperturbed model), Model 
B: H = 6 ■ 1CT 7 , Model C: X = 6 • 1CT 5 , Model D: H = 6 ■ 1CT 3 , 
Model E: H = 1.2 • 10~ 2 . A vertical E = Eo line is provided for a 
reference. 



suits derived without applying any scattering are indicated 
by N = 0. 



3 SPECTRA OF ACCELERATED PARTICLES 

In order to provide qualitative evaluations of turbulence ef- 
fects in the volume of reconnecting magnetic field, but con- 
sidering it only as a factor introducing random motion com- 
ponent to particle trajectories, we performed simulations of 
energetic proton spectra with a varying amount of turbu- 
lence (= scattering). As explained above, this approach as- 
sumes the existence of short wave magnetic field perturba- 
tions to be present in the limited volume - the black rect- 
angle in Fig. 1 - near the central neutral point, but we do 
not consider the influence of the turbulence on the reconnec- 
tion process. Thus it is a complementary approach to that 
using MHD modelling of the turbulent reconnection includ- 
ing wave perturbations from a narrow wave vector range, as 
discussed in section 1. 
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Figure 4. Acceleration times, T, versus the final particle energies, 
E, in the Model C (cf. Fig. 3). The simulation time cut-off at 
tmax = 100 s is indicated with a dashed line. Due to the applied 
trajectory splitting procedure the particles represented by points 
close to T = t max have very small weights and further increasing 
tmax does not lead to modification of the spectra presented in 
Fig. 3. 



We performed simulations of particle evolution starting 
at the same 'injection' energy Eo = 2 MeV, avoiding consid- 
eration of the real injection process at much lower energies 
(cf. Miller et al. 1997). For each set of particles we derived 
the spectrum of particles escaping from the reconnection 
volume, as illustrated in Fig. 3. In the non-perturbed model 
(N = 0, curve A) protons can increase their initial energy by 
approximately 70% . One should note that the injected en- 
ergetic particles can gain as well as loose energy. Introducing 
trajectory perturbations results in substantial modification 
of the acceleration process (curves B, C, D, E in Fig. 3). 
The spectrum energy cut-off shifts to higher values and the 
spectrum becomes harder when the amount of scattering 
('turbulence amplitude') is increased. In our simulations the 
resulting flat spectra extend up to 50 MeV for models with 
strong turbulence (Model D and E), and a steeper part of 
the spectrum is recorded at energies above 100 MeV. This 
behaviour results from the fact that the diffusive component 
introduced in to particle trajectories by the scattering en- 
ables some particles to stay in the reconnection region much 
longer and diffuse back close to the null point from outside. 
As illustrated, this can have a pronounced influence on the 
acceleration process by substantially increasing the parti- 
cle mean energy gain and providing much larger energies of 
individual particles. There is a general trend for the acceler- 
ation efficiency to increase with the perturbation amplitude 
in the considered range of values for the H parameter. 

Acceleration by the uniform electric field is character- 
ized with the energy-independent rate of particle energy in- 
crease, leading to a linear relation between the particle ac- 
celeration time and its final energy. In contrast, a plot of cor- 
relation between the individual particle acceleration time (= 
the physical simulation time) and its final energy presented 
in Fig. 4 shows that the acceleration rate depends on en- 
ergy, being preferentially defined by particle drifts in the 'V 
x B' electric fields. A large energy dispersion in a given time 
reflects variety of involved particle diffusive trajectories. 
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4 FINAL REMARKS 

Consideration of energetic particle acceleration accompany- 
ing the magnetic field reconnection process requires knowl- 
edge of the electromagnetic field structure in the region of 
interest. Until now the available numerical models of such 
fields were oversimplified and fail to consistently include to- 
gether the short and long MHD wave modes. Therefore, in 
the present simulations we apply an analytic model for a re- 
connecting field as the background, and the particle scatter- 
ing process is imposed as the small amplitude uncorrelated 
angular momentum perturbations. A comparison of the ac- 
celeration process in such a model to the unperturbed model 
reveals the important result that inclusion of the turbulent 
field component into the reconnection volume can change 
the acceleration process in a qualitative way, enabling par- 
ticles to reach much higher final energies and significantly 
increase their mean energy gain. 

A serious limitation on the validity of our simplified 
modelling arises due to a non-self-consistent introduction 
of the MHD turbulent field. The trajectory perturbations in 
the turbulent medium are considered without taking into ac- 
count the influence - in both directions - of these turbulent 
motions on the reconnection process itself. As mentioned in 
section 1, the presence of magnetic field turbulent structures 
increases mixing in the medium and is the source for anoma- 
lous resistivity leading to more effective reconnection as dis- 
cussed recently by Lazarian & Vishniac 2000. Subsequently, 
it leads to more efficient particle acceleration (Matthaeus 
et al. f984, Ambrosiano et al. 1988, Scholer & Jamitzky 
1989, Veltri at al. 1998, Birn & Hesse 1994, Kliem et al. 
1998, Schopper et al. 1999). As discussed in the present pa- 
per from a slightly different perspective it is expected to 
increase the cosmic ray acceleration efficiency, and influence 
the involved time scales and details of the resulting energy 
spectrum. However, for a given scattering conditions within 
the reconnection region the spectrum upper energy cut-off 
is limited by the 'global' perturbed structure of the recon- 
necting volume and not by the local conditions within the 
thin reconnection current sheet. 

The present work was supported by the Komitet Badan 
Naukowych through the grants PB 179/P03/96/11 and PB 
258/P03/99/17. 
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